Trabalho Final
Ana Flavia Arrym, Bruno Penãs, Isabela Brunele, Enzo Carmona, Gabriel Belisiario, Mariana Canto
2024-04-26
#Importando bibliotecas
library(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.0 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ ggplot2 3.4.1 ✔ tibble 3.2.0
## ✔ lubridate 1.9.2 ✔ tidyr 1.3.0
## ✔ purrr 1.0.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorslibrary(readxl)
library(dplyr)
library(clipr)## Welcome to clipr. See ?write_clip for advisories on writing to the clipboard in R.library(did)
library(ggplot2)
library(modelsummary)## `modelsummary` 2.0.0 now uses `tinytable` as its default table-drawing
## backend. Learn more at: https://vincentarelbundock.github.io/tinytable/
##
## Revert to `kableExtra` for one session:
##
## options(modelsummary_factory_default = 'kableExtra')
##
## Change the default backend persistently:
##
## config_modelsummary(factory_default = 'gt')
##
## Silence this message forever:
##
## config_modelsummary(startup_message = FALSE)setwd("C:/Users/Ana Flávia/OneDrive - Insper - Institudo de Ensino e Pesquisa/Economia/Quinto Semestre/Microeconomia IV/Trabalho Final")1 Ajustando as bases de dados
1.1 Base: mortes por overdose por estado
mortes = read.csv("NCHS_-_Drug_Poisoning_Mortality_by_State__United_States_20240417 (1).csv")
mortes = mortes %>%
select(State, Year, Sex, Age.Group, Race.and.Hispanic.Origin, Deaths, Population) %>%
mutate(Year = as.character(Year)) %>%
filter(Sex == "Both Sexes", Age.Group == "All Ages",
Race.and.Hispanic.Origin == "All Races-All Origins") %>%
select(State, Year, Deaths, Population) %>%
mutate(Overdoses = Deaths) %>%
select(State, Year, Overdoses, Population)%>%
rename(year = Year)
mortes <- mortes %>%
filter(State != "United States",
year != "1999",
year != "2000",
year != "2001",
year != "2002")
view(mortes)1.2 Base: IDH por estado
IDH = read.csv("GDL-Subnational-HDI-data.csv", sep = ";", quote = "",
header = TRUE)
IDH = IDH %>%
rename(State = Region) %>%
select(State, X2003, X2004, X2005, X2006,X2007, X2008, X2009,
X2010, X2011, X2012, X2013, X2014, X2015, X2016) %>%
slice(-1) %>%
pivot_longer(!State, names_to = "year", values_to = "IDH") %>%
mutate(year = gsub("X","", year))
view(IDH)1.3 Base: desemprego por estado
unemp = read_xls("emp-unemployment.xls", sheet = "States")## New names:
## • `` -> `...2`
## • `` -> `...3`
## • `` -> `...4`
## • `` -> `...5`
## • `` -> `...6`
## • `` -> `...7`
## • `` -> `...8`
## • `` -> `...9`
## • `` -> `...10`
## • `` -> `...11`
## • `` -> `...12`
## • `` -> `...13`
## • `` -> `...14`
## • `` -> `...15`
## • `` -> `...16`
## • `` -> `...17`
## • `` -> `...18`
## • `` -> `...19`
## • `` -> `...20`
## • `` -> `...21`
## • `` -> `...22`
## • `` -> `...23`
## • `` -> `...24`
## • `` -> `...25`
## • `` -> `...26`
## • `` -> `...27`
## • `` -> `...28`
## • `` -> `...29`
## • `` -> `...30`
## • `` -> `...31`
## • `` -> `...32`
## • `` -> `...33`
## • `` -> `...34`
## • `` -> `...35`
## • `` -> `...36`
## • `` -> `...37`
## • `` -> `...38`
## • `` -> `...39`
## • `` -> `...40`
## • `` -> `...41`unemp = unemp %>%
slice(-(1:4)) %>%
slice(-(54:62))
header = unemp %>%
slice(1) %>%
unlist()
unemp = unemp %>%
slice(-1) %>%
rename_with(~header) %>%
slice(-1) %>%
rename_with(
~ paste0("X", .)) %>%
rename(State = XArea,
Fips = XFips) %>%
select(State, X2003, X2004, X2005, X2006, X2007, X2008,
X2009, X2010, X2011, X2012, X2013, X2014, X2015, X2016) %>%
pivot_longer(!State, names_to = "year", values_to = "Desemprego") %>%
mutate(year = gsub("X","", year))
view(unemp)1.4 Base: renda por estado
income = read_xlsx("h08.xlsx")## New names:
## • `` -> `...2`
## • `` -> `...3`
## • `` -> `...4`
## • `` -> `...5`
## • `` -> `...6`
## • `` -> `...7`
## • `` -> `...8`
## • `` -> `...9`
## • `` -> `...10`
## • `` -> `...11`
## • `` -> `...12`
## • `` -> `...13`
## • `` -> `...14`
## • `` -> `...15`
## • `` -> `...16`
## • `` -> `...17`
## • `` -> `...18`
## • `` -> `...19`
## • `` -> `...20`
## • `` -> `...21`
## • `` -> `...22`
## • `` -> `...23`
## • `` -> `...24`
## • `` -> `...25`
## • `` -> `...26`
## • `` -> `...27`
## • `` -> `...28`
## • `` -> `...29`
## • `` -> `...30`
## • `` -> `...31`
## • `` -> `...32`
## • `` -> `...33`
## • `` -> `...34`
## • `` -> `...35`
## • `` -> `...36`
## • `` -> `...37`
## • `` -> `...38`
## • `` -> `...39`
## • `` -> `...40`income = income %>%
slice(-(1:60))
header = income %>%
slice(1) %>%
unlist()
income = income %>%
slice(-1) %>%
rename_with(~header) %>%
slice(-1) %>%
rename_with(
~ paste0("X", .)) %>%
rename(State = XState) %>%
select(State,X2003, X2004,X2005,X2006, X2007, X2008, X2009,
X2010, X2011, X2012, X2013, X2014, X2015, X2016) %>%
pivot_longer(!State, names_to = "year", values_to = "Income") %>%
mutate(year = gsub("X","", year))
view(income)1.5 Base: leis maconha por estado
leis = read_xlsx("BASE_LEISMACONHA.xlsx", sheet = "RML")
leis = leis %>%
select(-(2:32)) %>%
select(-(16:23))
leis = leis %>%
rename(State = "Estados/Anos") %>%
pivot_longer(!State, names_to = "year", values_to = "Tratamento") %>%
mutate(Tratamento = as.numeric(Tratamento))
view(leis)1.6 Base: uso de drogas (maconha e cocaína) por estado
drug_users = read_xlsx("Drug users.xlsx")
drug_users <- drug_users %>%
mutate(year = as.character(year))
view(drug_users)1.7 Base unificada
base_TF = left_join(mortes, IDH, by = c("State", "year"))
base_TF = left_join(base_TF, unemp, by = c("State", "year"))
base_TF = left_join(base_TF, income, by = c("State", "year"))
base_TF = left_join(base_TF, leis, by = c("State", "year"))
base_TF = left_join(base_TF, drug_users, by = c("State", "year"))
base_TF <- base_TF %>%
mutate(year = as.numeric(year),
weed_users = as.numeric(weed_users))
base_TF <- base_TF %>%
mutate(Fips = case_when(State == 'Alabama' ~ '01000',
State == 'Alaska' ~'02000',
State == 'Arizona' ~'04000',
State == 'Arkansas' ~'05000',
State == 'California' ~'06000',
State == 'Colorado' ~'08000',
State == 'Connecticut' ~'09000',
State == 'Delaware' ~'10000',
State == 'District of Columbia' ~'11000',
State == 'Florida' ~'12000',
State == 'Georgia' ~'13000',
State == 'Hawaii' ~'15000',
State == 'Idaho' ~'16000',
State == 'Illinois' ~'17000',
State == 'Indiana' ~'18000',
State == 'Iowa' ~'19000',
State == 'Kansas' ~'20000',
State == 'Kentucky' ~'21000',
State == 'Louisiana' ~'22000',
State == 'Maine' ~'23000',
State == 'Maryland' ~'24000',
State == 'Massachusetts' ~'25000',
State == 'Michigan' ~'26000',
State == 'Minnesota' ~'27000',
State == 'Mississippi' ~'28000',
State == 'Missouri' ~'29000',
State == 'Montana' ~'30000',
State == 'Nebraska' ~'31000',
State == 'Nevada' ~'32000',
State == 'New Hampshire' ~'33000',
State == 'New Jersey' ~'34000',
State == 'New Mexico' ~'35000',
State == 'New York' ~'36000',
State == 'North Carolina' ~'37000',
State == 'North Dakota' ~'38000',
State == 'Ohio' ~'39000',
State == 'Oklahoma' ~'40000',
State == 'Oregon' ~'41000',
State == 'Pennsylvania' ~'42000',
State == 'Rhode Island' ~'44000',
State == 'South Carolina' ~'45000',
State == 'South Dakota' ~'46000',
State == 'Tennessee' ~'47000',
State == 'Texas' ~'48000',
State == 'Utah' ~'49000',
State == 'Vermont' ~'50000',
State == 'Virginia' ~'51000',
State == 'Washington' ~'53000',
State == 'West Virginia' ~'54000',
State == 'Wisconsin' ~'55000',
State == 'Wyoming' ~'56000')) %>%
mutate(Fips = as.numeric(Fips))
base_TF<-base_TF %>%
mutate(taxa_overdose = (Overdoses/Population)*1000,
taxa_weed = (weed_users/Population)*1000,
taxa_cocain = (cocain_users/Population)*1000)
view(base_TF)2 Descritivas
2.1 Tabelas
#pré tratamento (2012)
base_TF %>%
filter(State %in% c("Colorado", "Nevada")) %>%
filter(year %in% c("2003":"2011")) %>%
rename(Cocaina = cocain_users,
Maconha = weed_users,
Renda = Income,
Habitantes = Population) %>%
select(State, Overdoses, Habitantes, IDH, Desemprego, Renda,
Cocaina, Maconha) %>%
datasummary(formula = ~ Overdoses + Habitantes + IDH
+ Desemprego + Renda + Cocaina + Maconha
~ State*(Mean + SD))| Colorado | Nevada | |||
|---|---|---|---|---|
| Mean | SD | Mean | SD | |
| Overdoses | 667.11 | 115.80 | 491.67 | 93.79 |
| Habitantes | 4807537.89 | 208834.96 | 2545901.44 | 169886.65 |
| IDH | 0.93 | 0.00 | 0.90 | 0.00 |
| Desemprego | 5.97 | 1.79 | 7.39 | 4.03 |
| Renda | 76971.11 | 3454.54 | 69088.89 | 4295.52 |
| Cocaina | 152024.49 | 19155.18 | 55942.56 | 7872.71 |
| Maconha | 713874.37 | 122214.88 | 302510.13 | 38672.94 |
#pré tratamento (2014)
base_TF %>%
filter(State %in% c("Alaska", "Oregon")) %>%
filter(year %in% c("2003":"2013")) %>%
rename(Cocaina = cocain_users,
Maconha = weed_users,
Renda = Income,
Habitantes = Population) %>%
select(State, Overdoses, Habitantes, IDH, Desemprego, Renda,
Cocaina, Maconha) %>%
datasummary(formula = ~ Overdoses + Habitantes + IDH
+ Desemprego + Renda + Cocaina + Maconha
~ State*(Mean + SD))| Alaska | Oregon | |||
|---|---|---|---|---|
| Mean | SD | Mean | SD | |
| Overdoses | 98.00 | 22.66 | 456.91 | 61.25 |
| Habitantes | 692375.27 | 29618.52 | 3748458.18 | 133901.36 |
| IDH | 0.93 | 0.01 | 0.91 | 0.01 |
| Desemprego | 7.19 | 0.54 | 7.88 | 2.04 |
| Renda | 79466.36 | 3991.66 | 65460.91 | 2490.47 |
| Cocaina | 16466.17 | 3102.24 | 83763.70 | 15332.42 |
| Maconha | 121482.87 | 18702.32 | 594795.02 | 100268.74 |
#summary geral
base_TF %>%
filter(!State %in% c("Alaska", "Oregon","Colorado", "Nevada")) %>%
rename(Cocaina = cocain_users,
Maconha = weed_users,
Renda = Income,
Habitantes = Population) %>%
select(State, Overdoses, Habitantes, IDH, Desemprego, Renda,
Cocaina, Maconha) %>%
datasummary(formula = ~ Overdoses + Habitantes + IDH
+ Desemprego + Renda + Cocaina + Maconha
~ Mean + SD)| Mean | SD | |
|---|---|---|
| Overdoses | 804.42 | 855.61 |
| Habitantes | 6280125.08 | 6965799.14 |
| IDH | 0.91 | 0.02 |
| Desemprego | 5.99 | 2.00 |
| Renda | 66166.26 | 10341.97 |
| Cocaina | 125089.51 | 151932.98 |
| Maconha | 699701.09 | 842269.89 |
2.2 Gráficos
#gráficos com os tratamentos
base_TF %>%
filter(State %in% c("Alaska", "Colorado", "Oregon",
"Washington")) %>%
filter(year %in% c("2003":"2013")) %>%
ggplot(aes(x = year, y = taxa_overdose))+
geom_line(color = "#336600")+
facet_wrap(~State) +
labs(title = "Overdose para os tratados",
x = "Ano", y = "Taxa de overdose por 1000 habitantes")
base_TF %>%
filter(State %in% c("Alaska", "Colorado", "Oregon",
"Washington")) %>%
filter(year %in% c("2003":"2013")) %>%
ggplot(aes(x = year, y = taxa_weed))+
geom_line(color = "#336600")+
facet_wrap(~State) +
labs(title = "Usuários de maconha nos grupos tratados",
x = "Ano", y = "Usuários de maconha por 1000 habitantes")
base_TF %>%
filter(State %in% c("Alaska", "Colorado", "Oregon",
"Washington")) %>%
filter(year %in% c("2003":"2013")) %>%
ggplot(aes(x = year, y = taxa_cocain))+
geom_line(color = "#336600")+
facet_wrap(~State) +
labs(title = "Usuários de cocaína nos grupos tratados",
x = "Ano", y = "Usuários de cocaína por 1000 habitantes")
#gráfico de geral
medias <- base_TF %>%
filter(!State %in% c("Alaska", "Oregon","Colorado", "Nevada")) %>%
group_by(year) %>%
summarise(medias_overdose = mean(taxa_overdose),
medias_maconha = mean(taxa_weed),
medias_cocaina = mean(taxa_cocain))
medias %>%
ggplot(aes(x = year, y = medias_overdose))+
geom_line(color = "#336600") +
labs(title = "Médias de overdose nos EUA - geral",
x = "Ano", y = "Médias de Overdose por 1000 habitantes")
medias %>%
ggplot(aes(x = year, y = medias_maconha))+
geom_line(color = "#336600") +
labs(title = "Médias de usuários de maconha nos EUA - geral",
x = "Ano", y = "Médias de usuários de maconha por 1000 habitantes")
medias %>%
ggplot(aes(x = year, y = medias_cocaina))+
geom_line(color = "#336600") +
labs(title = "Médias de usuários de cocaína nos EUA - geral",
x = "Ano", y = "Médias dos usuários de cocaína por 1000 habitantes")
3 Regressão
3.1 Regressão de morte por overdoeses
#tentando rodar o sttagerd
primeiros_tratamentos <- base_TF %>%
filter(Tratamento > 0) %>%
group_by(State) %>%
summarise(primeiros_tratamentos = min(year))
view(primeiros_tratamentos)
base_TF <- base_TF %>%
mutate(year = as.numeric(year),
Tratamento = as.numeric(Tratamento),
Fips = as.numeric(Fips)) %>%
left_join(primeiros_tratamentos, by = "State") %>%
mutate(Tratado = ifelse(Tratamento > 0, primeiros_tratamentos, 0)) %>%
arrange(State, year) %>%
fill(Tratado, .direction = "down") %>%
select(-primeiros_tratamentos)
base_TF <- base_TF %>%
mutate(Tratado = case_when(State == "Alaska"~"2014",
State == "California"~"2016",
State == "Colorado"~"2012",
State == "Maine"~"2016",
State == "Massachusetts"~"2016",
State == "Nevada"~"2016",
State == "Oregon"~"2014",
State == "Washington"~"2012",
.default = "0")) %>%
mutate(Tratado = as.numeric(Tratado))
view(base_TF)
#filtrei esses 4 estados pq não tem pós tratamento pra eles
base_2 <- base_TF %>%
filter(State != "California",
State != "Maine",
State != "Massachusetts",
State != "Nevada")
view(base_2)
overdose_death <- att_gt(yname = "taxa_overdose",
tname = "year",
idname = "Fips",
gname = "Tratado",
data = base_2)## Warning in pre_process_did(yname = yname, tname = tname, idname = idname, : Be aware that there are some small groups in your dataset.
## Check groups: 2012,2014.## Warning in att_gt(yname = "taxa_overdose", tname = "year", idname = "Fips", :
## Not returning pre-test Wald statistic due to singular covariance matrixaggte(overdose_death, type = "dynamic")##
## Call:
## aggte(MP = overdose_death, type = "dynamic")
##
## Reference: Callaway, Brantly and Pedro H.C. Sant'Anna. "Difference-in-Differences with Multiple Time Periods." Journal of Econometrics, Vol. 225, No. 2, pp. 200-230, 2021. <https://doi.org/10.1016/j.jeconom.2020.12.001>, <https://arxiv.org/abs/1803.09015>
##
##
## Overall summary of ATT's based on event-study/dynamic aggregation:
## ATT Std. Error [ 95% Conf. Int.]
## -0.0288 0.0055 -0.0395 -0.0181 *
##
##
## Dynamic Effects:
## Event time Estimate Std. Error [95% Simult. Conf. Band]
## -10 -0.0015 0.0025 -0.0077 0.0048
## -9 -0.0092 0.0067 -0.0260 0.0075
## -8 -0.0004 0.0086 -0.0221 0.0212
## -7 -0.0002 0.0052 -0.0131 0.0128
## -6 0.0080 0.0253 -0.0554 0.0714
## -5 0.0104 0.0038 0.0008 0.0201 *
## -4 -0.0227 0.0222 -0.0784 0.0331
## -3 0.0050 0.0066 -0.0115 0.0215
## -2 -0.0093 0.0145 -0.0457 0.0272
## -1 -0.0096 0.0122 -0.0403 0.0211
## 0 0.0025 0.0049 -0.0099 0.0149
## 1 -0.0126 0.0037 -0.0219 -0.0033 *
## 2 -0.0327 0.0098 -0.0573 -0.0081 *
## 3 -0.0351 0.0074 -0.0538 -0.0164 *
## 4 -0.0660 0.0101 -0.0914 -0.0405 *
## ---
## Signif. codes: `*' confidence band does not cover 0
##
## Control Group: Never Treated, Anticipation Periods: 0
## Estimation Method: Doubly Robust#gráficos
ggdid(aggte(overdose_death, type = "dynamic"))
3.2 Regressão de uso de maconha
weed <- att_gt(yname = "taxa_weed",
tname = "year",
idname = "Fips",
gname = "Tratado",
data = base_2)## Warning in pre_process_did(yname = yname, tname = tname, idname = idname, : Be aware that there are some small groups in your dataset.
## Check groups: 2012,2014.## Warning in att_gt(yname = "taxa_weed", tname = "year", idname = "Fips", : Not
## returning pre-test Wald statistic due to singular covariance matrixaggte(weed, type = "dynamic")##
## Call:
## aggte(MP = weed, type = "dynamic")
##
## Reference: Callaway, Brantly and Pedro H.C. Sant'Anna. "Difference-in-Differences with Multiple Time Periods." Journal of Econometrics, Vol. 225, No. 2, pp. 200-230, 2021. <https://doi.org/10.1016/j.jeconom.2020.12.001>, <https://arxiv.org/abs/1803.09015>
##
##
## Overall summary of ATT's based on event-study/dynamic aggregation:
## ATT Std. Error [ 95% Conf. Int.]
## 6.493 6.7584 -6.7533 19.7393
##
##
## Dynamic Effects:
## Event time Estimate Std. Error [95% Simult. Conf. Band]
## -10 -4.7907 2.8072 -12.9850 3.4037
## -9 10.9651 3.6991 0.1675 21.7628 *
## -8 -11.8953 4.8048 -25.9207 2.1300
## -7 -5.3488 5.7925 -22.2571 11.5595
## -6 9.2674 1.4211 5.1194 13.4155 *
## -5 4.8605 3.5820 -5.5955 15.3164
## -4 8.0872 1.3330 4.1962 11.9783 *
## -3 0.7674 3.3248 -8.9377 10.4726
## -2 13.7616 3.4050 3.8225 23.7007 *
## -1 1.9593 8.1855 -21.9343 25.8529
## 0 -8.6802 3.9383 -20.1763 2.8158
## 1 -1.8370 14.1333 -43.0920 39.4181
## 2 9.0236 9.7359 -19.3955 37.4427
## 3 14.6732 14.2983 -27.0636 56.4099
## 4 19.2854 9.9434 -9.7395 48.3103
## ---
## Signif. codes: `*' confidence band does not cover 0
##
## Control Group: Never Treated, Anticipation Periods: 0
## Estimation Method: Doubly Robust#gráficos
ggdid(aggte(weed, type = "dynamic"))
3.3 Regressão de uso de cocaína
cocaine <- att_gt(yname = "taxa_cocain",
tname = "year",
idname = "Fips",
gname = "Tratado",
data = base_2)## Warning in pre_process_did(yname = yname, tname = tname, idname = idname, : Be aware that there are some small groups in your dataset.
## Check groups: 2012,2014.## Warning in att_gt(yname = "taxa_cocain", tname = "year", idname = "Fips", : Not
## returning pre-test Wald statistic due to singular covariance matrixaggte(cocaine, type = "dynamic")##
## Call:
## aggte(MP = cocaine, type = "dynamic")
##
## Reference: Callaway, Brantly and Pedro H.C. Sant'Anna. "Difference-in-Differences with Multiple Time Periods." Journal of Econometrics, Vol. 225, No. 2, pp. 200-230, 2021. <https://doi.org/10.1016/j.jeconom.2020.12.001>, <https://arxiv.org/abs/1803.09015>
##
##
## Overall summary of ATT's based on event-study/dynamic aggregation:
## ATT Std. Error [ 95% Conf. Int.]
## 0.0943 0.9362 -1.7406 1.9292
##
##
## Dynamic Effects:
## Event time Estimate Std. Error [95% Simult. Conf. Band]
## -10 0.1968 1.1379 -2.4435 2.8371
## -9 0.2638 0.8841 -1.7875 2.3151
## -8 -1.3440 1.4136 -4.6239 1.9359
## -7 0.9081 2.6160 -5.1617 6.9780
## -6 3.7719 2.5132 -2.0594 9.6033
## -5 0.4120 1.5634 -3.2155 4.0394
## -4 -0.4363 2.8796 -7.1176 6.2450
## -3 -0.2460 1.4337 -3.5726 3.0806
## -2 -0.3545 1.1041 -2.9164 2.2074
## -1 -1.6069 0.7456 -3.3368 0.1230
## 0 -0.5664 1.2559 -3.4805 2.3477
## 1 1.2170 1.5075 -2.2808 4.7149
## 2 2.6266 1.3369 -0.4754 5.7286
## 3 -0.3460 1.4724 -3.7623 3.0704
## 4 -2.4598 0.4880 -3.5920 -1.3276 *
## ---
## Signif. codes: `*' confidence band does not cover 0
##
## Control Group: Never Treated, Anticipation Periods: 0
## Estimation Method: Doubly Robust#gráficos
ggdid(aggte(cocaine, type = "dynamic"))